function S = cicaarsep(invA0,Atau,X,intercept)
% CICAARSEP separate sources using CICAAR parameters
%
%   S = cicaarsep(invA0,Atau,X [,intercept] )
%   S = cicaarsep([],A,X [,intercept] )
%
% note: intercept is optional
%
% Copyright: Mads Dyrholm, 2007
if isempty(invA0)
  invA0 = pinv(Atau(:,:,1));
  Atau = Atau(:,:,2:end);
end
K = size(invA0,1);
[D,N] = size(X);
if nargin>3
  X = X - repmat(intercept,[1,N]);
end
if isempty(Atau)
  S = invA0*X;
  return
end
[D,KC,L] = size(Atau);
S = zeros(K,L+N);
AtauMAT = reshape(Atau,[D KC*L]);
invA0AtauMAT = invA0*AtauMAT;
invA0X = invA0*X; 
for t=1:N
  sVEC = S(1:KC,t-1+L:-1:t);
  S(:,t+L) = invA0X(:,t) - invA0AtauMAT * sVEC(:);
end
S = S(:,L+1:end);
